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ABSTRACT 

^  This  paper  is  concerned  with  one-step  dif ference  methods  for  parabolic 
history  value  problems  in  one  spaoe  variable.  These  problems,  which  have  the 
feature  that  the  evolution  of  the  solution  is  influenced  by  'all  its  past' 
occur  in  the  theory  of  viscoelastic  liquids  (materials  with  'memory').  The 
history  dependence  is  represented  by  a  Volterra-integral  in  the  equation  of 
motion.  Using  recently  obtained  existence  results  <see  Henardy  (1981b))  and 
smoothness  assumptions  on  the  solution,  we  derive  a  local  stability  and 
convergence  result  for  a  Crank-Ni col son- type  difference  scheme  by  interpreting 
the  linearized  scheme  as  perturbation  of  a  strictly  parabolic  scheme  without 
memory  term.  Second  order  convergence  is  shown  on  sufficiently  small  time 
intervals.  The  presented  approach  carries  over  to  other  one-step  difference 
methods  like  implicit  and  explicit  Euler  schemes. 
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SIGNIFICANCE  AND  EXPLANATION 


In  this  paper  we  deal  with  the  numerical  solution  of  parabolic  history 
value  problems*  These  problems  have  the  feature  that  the  governing  equation 
depends  on  the  history  of  the  solution  such  that  it  is  posed  as  functional 
differential  equation  (that  means  that  the  equation  can  involve  Volterra  type 
integrals  and  not  only  derivatives  of  the  function  in  question).  Problems  of 
this  kind  occur  in  the  theory  of  viscoelastic  fluids  and  there  the  functional 
term  of  the  equation  represents  the  'memory*  of  the  material.  We  devise  a 
finite  difference  method  for  the  numerical  solution  of  such  problems  and 
investigate  the  convergence  properties.  It  turns  out  that  this  method  (which 
is  of  Crank-Nicolson  type)  is  second  order  accurate  as  the  grid  parameters 
tend  to  zero. 


The  responsibility  for  the  wording  and  views  expressed  in  this  descriptive 
summary  lies  with  MRC#  and  not  with  the  author  of  this  report. 


DIFFERENCE  METHODS  FOR  PARABOLIC  HISTORY  VALUE  PROBLEMS 

* 

Peter  Markowich 

1 .  Introduction 

We  are  concerned  with  finite  difference  methods  for  scalar  functional  differential 
equations  with  prescribed  history  data: 

(1.1)  utt  =  g(uxt,u^),  x  e  1-1.1),  t  >  0 

(1.2)  u(x,t>  *  u(x.t),  x  e  I-1.1J,  t  *  0 

(1.3)  g(uxt(-l,t),  ux(-1,*))  =  f_(t),  t  >  0 

(1.4)  g(uxt(1,t),  ux(1,*))  =  f+(t),  t  >  0  . 

Here  we  denoted  the  history  of  a  function  y  e  C(t-~#0])  (which  is  the  space  of  all 

continuous  functions  on  (-*,0]  with  a  finite  limit  at  t  =  -•)  by  yt(s)  =  y(t+s), 

s  4  0  and  the  (possibly)  nonlinear  functional  g  :  R  *  ft  ♦  R  where  SI  is  an  open  set  in 

C( [-*,0] ).  We  also  assume  that  the  Frechet  derivative  (with  respect  to  the  first  argument) 

O^g  >  €  >  0  in  R  *  ft.  This  assumption  makes  it  possible  to  interpret  (1.1)- (1.4)  as  a 

parabolic  initial-value  problem  in  a  certain  Banach  space  (see  Renardy  (1981b)).  Histories 

for  the  boundary  data  u  (-1,t)  and  u  (1,t)  are  u  (-1,t)  and  u  (1,t)  resp. 

f*  —  x  x 

Assuming  that  (1.3),  (1.4)  with  the  corresponding  prescribed  histories  can  be  solved 
for  ux(-1,t)  and  ux(1,t)  resp.  we  are  left  with  a  (parabolic)  history  value  problem 
with  Neumann  boundary  conditions  at  x  =  ±1. 

Problems  of  this  kind  occur  in  the  theory  of  viscoelastic  liquids  when  the 
constitutive  law  is  expressed  as  a  function  of  the  strain  history  (see  Lodge  (1974),  Lodge, 
McLeod  and  Nohel  (1978)  and  Renardy  (1981b)).  The  functional  g  is  assumed  to  be  of 
Volterra  type: 

<1.5)  g(uxt'ux>  ”  Y(uxt'Ux>  +  /-»  a(t-s>b<V(t)'ux(s>,ds  * 

2 

where  Y,  b  :  R  ♦  R  and  a(<J)  is  an  exponentially  decaying  (as  0  ♦  •)  memory  kernel. 
The  dependence  of  u  on  the  space  variable  x  is  not  stressed  explicity. 

•  ■  —  1 
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In  this  paper  we  set  up  a  Crank-Nicolson-type  discretisation  for  (1.1)  on  an 
interval  (0,T)  with  given  history  (1.2).  The  boundary  problem*  (1.3),  (1.4)  are 
discretized  by  the  mid-point  rule  thus  giving  discrete  Neumann  boundary  data  for  the  Crank - 
Nicolson  scheme.  Assuming  that  T  is  sufficiently  small  and  that  the  solution  u  of 

(1.1 ) -(l .4)  is  sufficiently  smooth  we  show  stability  of  the  linearized  difference  scheme 
and  consistency  of  the  nonlinear  scheme  at  the  exact  solution  u  in  a  discrete  Sobolev 
space  norm.  From  this  and  from  the  uniform  (in  the  mesh-sizes)  Lipschitz  continuity  of  the 
linearized  scheme  we  conclude  convergence  (of  order  two)  from  Keller's  (197$)  theory. 

The  approach  is  to  interpret  the  linearized  difference  scheme  as  perturbation  of  a 
strictly  parabolic  scheme  (without  history  term)  and  stability  of  the  scheme  for  the 
history  value  problem  will  be  concluded  from  the  stability  of  the  parabolic  difference 
scheme.  Therefore  this  approach  is  applicable  to  other  one-step  difference  schemes  like 
the  implicit  and  explicit  Euler  schemes.  The  implicit  Euler  scheme  may  be  chosen  if 
approximations  are  needed  on  a  large  interval  [0,T]  (assuming  the  exact  solution  exists 
there)  because  it  is  strongly  A-e-stable  (see  Markowich  and  Renardy  (1981a, b)). 

The  paper  is  organized  as  follows.  In  Section  two  we  define  the  function  spaces  which 
w©  will  need  and  introduce  some  notations.  Section  three  deals  with  the  discretization  of 
the  boundary  problem  and  Section  four  is  concerned  with  the  Crank-Nicolson  scheme  for 
(1.1  ). 

2.  Definitions  and  Notations 
We  denote 

(2.1)  Ci([-»#t1l)  -Ci((-*\tinn  ♦  R|li»  f(t)  is  finite  and 

lim  f  ^  * ( t )  «  0  for  j»1,...,i),  ieM  • 

*►-»  0 

for  some  t-,  8  R«  For  y  Q  we  define  the  history  yt  8  CA((-*,Oj)  by 

yt(s)  -  y ( t+s ) ,  s  8 
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Renardy  (1981a)  investigated  a  model  for  the  stretching  of  a  filament  of  polymeric 


liquid  in  the  form  of  an  initial  value  problem  ((1*1),  (1.3),  (1.4)  are  assumed  to  hold  on 
T),  T  6  R  U  {»}  and  u(x,t  *  -•)  is  prescribed).  He  used  a  functional  g  of  the 
form  (1.5).  A  global  (T  =  •)  existence  and  uniqueness  theorem  for  sufficiently  small 
data  f+,  f_  (in  the  sense  of  a  certain  Sobolev  space)  and  a  local  (for  T  <  0,  |t| 
large)  existence  and  uniqueness  theorem  for  arbitrarily  large  boundary  data  was 
established*  The  boundary  problem  (1.3)  (in  initial  value  form)  was  investigated 
analytically  and  numerically  in  Markowich  and  Renardy  (1981a)  using  an  implicit  Euler-type 
discretization. 

The  full  spatial-temporal  problem  (in  initial  value  form)  was  investigated  numerically 
in  Markowich  and  Renardy  (1981b).  Again,  an  implicit  Euler-type  discretization  was  usc^  in 
order  to  get  approximate  solutions  with  the  same  asymptotic  behaviour  as  the  exact  solution 
(as  t  ♦  •)  and  in  order  to  cope  with  the  singular  perturbation  character  of  the  problem 
(the  Newtonian  contribution  of  the  viscosity  acts  as  singular  perturbation  parameter). 

Lodge,  McLeod  and  Nohel  (1978)  investigated  the  history  value  problem  for  the  boundary 
problem  (1.3)  assuming  the  relation  (1.5)  with  Y  specified.  Under  certain  assumptions 
on  a  and  b  they  proved  a  global  existence  theorem  and  investigated  the  asymptotic 
behaviour  of  the  solution  as  t  ♦ 

Nevanlinna  (1978)  employed  an  implicit  Euler-type  discretization  for  the  boundary 
history  value  problem  and  proved  uniform  convergence  on  of  the  order  O(h^), 

0  <  \  <  1. 

Renardy  (1981b)  proved  a  local  (for  t  6  10, T],  T  sufficiently  small)  existence  and 
uniqueness  theorem  for  the  history  value  problem  (1.1)-(1.4)  under  mild  assumptions  on  the 
functional  g  (not  using  the  special  form  (1.5)).  He  transformed  (1.1)  to  a  system  of 
equations  which  can  be  interpreted  as  parabolic  in  the  Banach  space  C([-*,0),  ( (-1 , 1 1) )• 

This  system  can  be  treated  following  Sobolevskii  (1966), 
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In  the  sequel  g  always  denotes  a  functional  from  R  *  ft  into  R  where  ft  Is  an 
open  set  in  some  space  of  functions  which  are  defined  on  [*-*#0].  Frechet  derivatives  of 


g  are  denoted  by  indices  ( denotes  the  derivative  of  g  with  respect  to  the  first 
argument,  gj  with  respect  to  the  second  argument,  for  example  denotes  the  second 

Frechet  derivative  of  g  obtained  by  differentiating  first  with  respect  to  the  first 

argument  and  then  with  respect  to  the  second  argument)* 

2  £ 

L  <[-1,1])  denotes  the  space  of  square  integrable  functions  on  [-1,1],  H  ([-1,1]) 

denotes  the  space  of  L^U-1,1])  functions  whose  (generalized)  derivatives  of  order  up  to 

1  are  square  integrable  and  C™( [t  t  ],  H*( 1-1,1])  denotes  the  space  of  Cm-functions 

1  ,  z 

u  i  ♦  H*([-1,1]).  All  these  spaces  are  equipped  with  their  natural  norms. 

For  the  difference  scheme  we  define  a  grid  Gr(h,k)  on  the  infinite  strip 
[-1,1]  *  [-*, T]  by  setting  h  «  k  *  J;  H,  N  e  M  such  that 


(2.2) 


Gr  (  h ,  k  )  -  {(x  ,t  ) | x , *  jh,  j  -  (-M-1 ) ( 1 ) (M+1 ) ;  t  *  nk,  n  <  N) 
j  n  j  n 


holds.  The  exterior  gridpoints  ( x_M- 1 , tn ) ,  (x^jt .  )  will  be  needed  to  define  a  second 

order  approximation  to  the  Neumann  boundary  conditions.  We  denote  grid  functions  by 

U(N)  *  ((U*)N  _)*+*  0^  e  R.  We  also  need  grid  functions  on  Gr(h,k)  -  {exterior  grid 

3  n*  j— w-i  i 

points) «  They  are  defined  analogously. 
f(N) 


By  we  denote  the  grid  function  on  Gr(h,k)  O  {(x^,tn)| j  *  i,  n  <  n): 


(2.3) 


u.(N)  -  <u")N 

i  in* 


(N) 


For  the  discretization  of  the  boundary  problems  we  need  the  {t  *  nkjn  <  N)  grid  on  the 

n 

real  line  and 
nd 

(2.4)  y 

are  functions  on  this  grid.  Also  we  set  We  define  the  following  time- 

difference  quotients: 

(2.5 )  (a)  o  y 


/  n.N  n  _  _ 

<y  ’n— •  /ye* 


,  n*M  n 

sV  -  * — 


(2.5) (b) 


n  n-1 

6’yn  .  1^11 - 

y  k 


i- 

i 

i 

i 

i 


{ 2  *5 ) { c) 


5y 


2k 


Obviously  6+y°  =  6“yn+1,  6yn  *  ^(6+yn  +  bolds.  We  need  the  spatial  differences? 


( 2.6) ( a) 


Vi 


(2.6) ( b) 


Vi 


(2.6) (c) 


Az. 


i  2h 

These  difference  quotients  will  also  be  applied  (component  wise)  to  grid  functions,  for 

i  Ai1(N)  /Ai!niN 

example  Aik  •  AUi  n*-°“* 


3.  The  Boundary  Problem 

In  order  to  solve  the  boundary  problems  (1.3),  (1.4)  for  ux(1,t)  and  ux(-1,t) 
resp.  we  discretize 

(3.D  g(y' (tj.y1)  -  f(t)  ,  t  e  [o,t]  ,  Ten 

(3.2)  y(t)  -  y(t)  ,  t  e  [-*,0] 
by  the  midpoint  rule 

(3.3)  g(«V\  (iky(N))tn+V2)  -  f(t„+1/2)  =  0,  0  <  a  <  M 


(3.4) 


y  -  y(tn)  -  0  ,  n  <  0 

(N)  ,  , 


Here  yn  denotes  the  approximation  to  y(tn),  y'  *  (y  )naB_„  and  i^  Is  the  linear 


interpolation  operator  defined  by 

.(N)  _  ,..nxN  e  i  # 


( 3.5 )  ( a)  i  ?  *  {u  NJ  =  <u  ) e  *'  lim  u  *  u  e  R>  — ►  C( [-»,T]  ) 

*  N  ri=”  i  >.00 


( 3.5 ) (b)  (iku^N^)(t)  »  u*  +  6+ui*(t-ti)  for  ti  4  t  <  t^+1 ,  i  <  N  . 

We  denote  the  right  hand  side  of  (3.3),  (3.4)  by 

(N) 


(3.6) 


Fk(y  ')  -  0 


where 


(3-7>  Fk  *  flN  *  ®N  '  SN  °Pen'  aN  C  \  * 

Here  A^  is  equipped  with  the  norm 

(3.8)  lu(N)l  *  maxjun|  +  max|6  un|  +  max  1 6*6  u°| 

N  n<N  n<N  1<n<N 

B^j  equals  A^  (as  a  set)  but  as  norm  we  take 

(3.9)  *w(N).  *  max|wn|  +  max|6  wn| 

N  n<N  n<N 

ft^  is  open  in  AN  since  g(v,»),  v  e  R  is  defined  on  an  open  set  ft  c*  C(  [-•,()]). 

The  convergence  analysis  of  (3.6)  proceeds  along  the  lines  of  Keller's  (1975) 
stability-consistency  concept.  For  consistency  we  need  smoothness  assumptions  on  (3.1), 
(3.2).  We  assume  that  there  is  locally  unique  solutions  y  of  (3.1),  (3.2)  and  that 

(3.10) (a)  g  e  C2(m  x  0,  R),  a  c  c(  [--,0]  ),  g1  >  e  >  0  on  R*0 

(3.10)  (b)  fec3([0,T)) 

(3.10)  (c)  y  e  c4([0,T]),  y  e  C4([-»,0]); 

(3. 10 )  ( d)  y  e  C2([-»,T)),  y*  6  8  for  t  e  [0,T] 

holds.  A  local  existence  (and  uniqueness)  theorem  is  given  in  Renardy  (1961b).  (3.10d) 

holds  if  the  compatibility  requirements 

(3. 1 1 )  (a)  g(y'(0),  y)  -  f(0) 

(3. 1 1 )  (b)  g1 (y1 (0) ,  y)J-(0)  +  g2(y'(0),  y)y •  -  f»(0) 
are  fulfilled. 

We  check  consistency  of  the  scheme  (3.3),  (3.3).  We  denote 

(3.12)  y<N)  -  <y(t  )>"  .  . 

n  n*-" 

( 3. 10 )  (c ) ,  (d)  imply  that  (i^y^)  e  ft  for  k  sufficiently  small  and  so  6 

holds.  Obviously  (F^(y^N^))n  »  6  (F^(y^N^))n  »  0  for  n  <  0.  Since 

(3.13)  -  ylt-,T)  "  0(lt2) 

-  2  2 

(because  y  e  C  ((-*,0J),  y  6  C  ((0,Tj)  and  t  *  0  is  a  grid  point)  we  obtain  using 

(3.10) (a),(b)  and  the  mean  value  theorem 

(3.14)  l(Fk(y<N>))n|  -  0(k2)  . 

The  smoothness  assumptions  (3.10)(c),(d)  imply  that 
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(3.15) 


( .  {  j(.)  Vl/2 

■  k  ,  .  v  n.  ^ ,,  2  „ 

' - k - <*'*  '[— ,0]  =  0(k  > 


holds.  Using  (3.10),  (3.15)  and  the  differentiated  equation  (3.1) 

(3.16)  g,|  (y'  (t)  ,yfc)y"(  t)  +  g2(y' (t) ,yfc) (y* ) -  f'(t)  =  0 
((3.10)(d)  implies  that  (y')t  e  C*([-*,01))  we  get,  after  a  simple  calculation 

(3.17) 

Since  all  estimates  hold  uniformly  for  n  <  N  we  obtain  from  (3.14),  (3.17) 

(3.18) 


|«+<Fk(ytN)))n|  »  0(k2)  ,  n  <  N 


IF.  (y(N))l  «  o(k2) 

k  bn 


Therefore  the  scheme  (3.3),  (3.4)  is  consistent  of  order  2  at  y'N^  if  (3.10)  holds. 


A  (  N  ) 

To  check  stability  we  calculate  the  Frechet-derivative  of  F^  at  y  and  get  for 


u(N)  e  ^ 


\  ,„t,A(N>,  (NKn+1  a(N).  n*H/2,e+  n 

(3.19(a)  (F£(y  ) u  )  =*  g1  ( o  ytt^Mr^y  )  )6  u  + 


.  ,*+  *  ( N ) . fcn+ 1 /2 w  ,  (N ) . fcn+1/2  .  £ 

+  g2(°  yCt^Mi^y  )  Mi^u  '  0  <  n  <  N 


(3.19) (b) 


,a(N)  .  (N ) .  n  n 
(F^(y  )  u '  )  =*  u  ,  n  <  0  . 


In  order  to  get  an  estimate  on  the  norm  of  (F/(y^))  j  B„  ♦  A.  we  investigate  the 

k  N  N 


equation 

(3.20) 


,  . a (N )  .  (N)  (N )  (N)  ,,-i.N 

F^(y  ) u  3  4  9  *  *  (*  ®  Bn  . 


For  n  <  0  we  have  un  *  v?n,  6  u11  =*  6  and  for  0  <  n  <  N 


(3.21) 


.+  n  a  (N) .  n+1/2  n+1 

0  u  =  gR(k) ( i^u  )  +  ^ 


where  gn(k)  is  the  linear  functional 


(3.22) 


gn(k) 


n+1/2> 

g1(^y(tn),(ik$(N))tn",/2) 


and 
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, .+*  .  .  .  .  •'(N)  n+1/2 \ 

g. ($  y(t  >, (1  y  )  ) 

I  n  K 


holds.  The  assumption  g^  >  e  >  0  on  R*fl  implies  that  gn(k)  is  uniformly  bounded 


(in  k,n)  and  that 


|5 n+  |  <  const|^n  |  uniformly  in  k,  n  <  N 


holds*  Prom  (3.21)  we  get: 


|un+^|  <  |un|  +  ck  max  |u^|  +  ck  | 

1 <j<n+l 


and  setting  vn  *  max  |un|,  2n  »  max  we  obtain 

1<j<n  j<n 


n+1  ^  n  .  n+1  t  n+1 

v  <  v  +  ckv  +  ckz 


n+1  .In  ,  n+1 
v  <  - — -  v  +  c  kz 
1-ck  1 


immediately  get  |vn|  <  c  ecTT  max  |z^|  and  therefore 

1<j<n 


holds.  We  obtain  from  (3.21) 


max  |un|  <  const,  max  ]^n| 


max  |6+un|  <  const,  max  |^n| 


In  order  to  prove  stability  in  the  norms  of  we  apply  $  to  (3.21)  and  obtain 

#i«,W  .  a  w,Vv2  5  (k).3  (K,  ...  t  1/5 

(3.29)  «“«V  -  gn(k)(— - jj - - -  )  +  (— - - )(iku(N))  n'1/2 

+  6*Vn  ,  n  >  0  • 

Vk)~Vi(k) 

(3.10)  implies  that  - j -  is  a  uniformly  bounded  functional  (in  k,n)  and 

therefore  the  second  term  on  the  right  hand  side  of  (3.29)  is  bounded  by  const,  max  |^n( 

n<N 

We  get  by  a  simple  calculation  (similar  to  (3.15)) 


UOjW  .  u  AH)  , "n-l /2 


(3.30)  I- 


10  -1  *  max  l6+u,1|  +  \  “ax  l4"«fun| 

l0*  1  n<N  i  n<N 


t 


Repeating  the  argument  following  (3.21)  we  get 
(3.31 ) 


max  |<5  6  u  |  4  const.  ( max  1 6  $  |  +  raax|^  |) 
1<n<N  n<N  n<N 


and  stability  follows: 
(3.32) 


(F^(y<N) ) )  1«B  <  const. 


where  const,  is  independent  of  k.  We  also  have  to  show  uniform  Lipschitz  continuity  of 
F^  in  a  neighborhood  S^(y^N^)  c  For  a^N^#  b^“ '  e  sp(y^*)  with  P  sufficiently 


small  we  get  from  (3.19)  and  the  mean  value  theorem: 
(3.33) 


ir(a(N))  -  FMb(N)),u(N).  <  const..a(N) 

k  k  BN 


(N)  . 


\« 


.  (N)i 

lu  l 


for  all  u(N)  e  where  const,  is  independent  of  k.  (3. 10) (a)  was  used  for  (3.33). 
Therefore  uniform  Lipschitz  continuity  of  F^  holds  on  Sp(y^N^). 

Now#  having  proven  consistency#  stability  and  Lipschitz  continuity  we  apply  Keller’s 
(1975)  theory  which  gives: 

Theorem  3.1.  Let  the  assumptions  (3.10)  on  the  history  value  problem  (3.1),  (3.2) 

hold.  Then  for  all  k  sufficiently  small  there  is  a  locally  unique 
(N) 


(3.34) 


solution  y  e  fl  Avt  of  the  midpoint  rule  (3.3),  (3.4)  and 
N  N 

max|yn-y(t  )|  +  max|6yn  1  -  y'(t  )|  + 
n<N  n  n<N  n_1 


+  max|6+6  yn  '  -  y"(t  )|  =  O(k^)  as  k  ♦  0 
n<N  n_1 


Moreover  the  (abstract)  Newton  procedure 


(3.35)  (a) 

( 3 . 35 )  ( b) 


(N)  =  -(H) 

(O)  (0) 


(N) 


(N) 


( i+1 ) 


(i) 


(F’(v(N)))-1F  (v(N)) 

1  k1  <i)M  V  (i)  ' 


i  >  0 

converges  quadratically  (to  y^)  from  a  sphere  of  starting  values 
which  does  not  shrink  as  k  ♦  0. 

The  abstract  Newton  method  (3.35)  immediately  translates  into  the  Newton  method  for 
determining  y\...#yN  from  (3.3)  assuming  that  (y^)^  m  is  given. 
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For  the  convergence  analysis  of  the  boundary  problem  we  do  not  require  T  small 
(except  that  a  smooth  solution  of  (3.1  i.  (3.2)  has  to  exist  on  l-*,T]).  Therefore  Theorem 
(3.1)  holds  for  any  finite  T  to  which  the  solution  y  can  be  smoothly  continued. 

However,  the  stability  constant  (3.32)  depends  on  T  and  so  does  the  estimate  (3.34), 

4.  The  Parabolic  Problem 


we  now  discretize  the  full  spatial-temporal  problem: 


(4.1) 

utt  -  «<uxt'ui,x  ■  °- x  ®  1  ® 

(0,Tl 

(4.2) 

u(x,t)  -  u(  X ,  t )  ,  X  €  1-1,1],  t  e  I-oo, 

0] 

(4.3) 

uxd  ,t)  -  y+(t),  t  s  t--,T) 

(4.4) 

ux(-i,t)  «  y_(t),  t  e  I-*,t) 

where  y+,  y_ 

solve  the  boundary  problems 

(4.5) 

g(y^(t),  y*)  «  f±(t),  t  a  [o,t] 

(4.6) 

y±(t)  -  ux(±i,t),  t  e  l--,o) 

which  fulfill  the  assumptions  of  Section  3.  For  convenience  we  carry  out  the 
differentiation  in  (4.1)  (assuming  sufficient  smoothness) 

(4.7)  utt  -  <g,<uxt,uXxt  +  92(uxt'ux)uxx>  “  °'  X  8  M',)'  1  8  t0*T1  • 

We  discretize  (4.7),  (4.2),  (4.3),  (4.4)  by  the  Crank-Nicolson  method: 

(4.8)  a+«"(^  -  (g^Ait^,  (ikAU*N))tn)«A+A_u£  ♦ 

♦  g2(Afiu£,  (ikAU*N))  n)(ikA+A_u|N>)  ")  -  0 
for  1  -  0  <  n  <  M 

(4.9)  u"  -  u(x.  ,t  )  -  0,  i  -  -M(1)M,  n  <  0 

i  in 

(4.10)  Aujj  ‘  y+  “  0'  n  *  N 

(4.11)  AU^  -  y"  »  0,  n  <  N  . 

Uj  denotes  the  approximation  to  u(xi#tn)  end  *  (U')^^.  The  boundary  values 

y^,  n  <  N  are  computed  by  discretizing  (4.5),  (4.6)  according  to  Section  3  and 
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therefore  are  assumed  to  be  known.  In  (4.10),  (4,11)  we  introduced  the  exterior  grid 


(4.18) 


cSx)  -  {z*N)  »  (Zj)N--(JZj  e  X,  lim  e  x> 

*  3“~  j+-» 

lz(N)l  -  max  lznl  ♦  max  l6“zn»  . 

c’(X)  n*N  n<N 


We  regard  G,  .  as  the  following  mappings 

K  /  n 


(4.19) 


G  .  :  «  e  C*(H*)  ♦  C  (H^#  L*  )  ( 

k,h  N  ,M  k  h  k  h  h 


assuming  that  f^(t)  fulfill  (3.10) (b)  and  that  the  boundary  values  y^(t)  fulfill 
( 3. 10 ) (c ) , (d) .  Moreover  g  e  C3 (R  *  ft,  R),  ft  c  C(  [-•,()}  )  open, 
g  e  C3(B  x  a,  L2<  (-1,1] )),  a  c  c(  [--,0] ,  L2<[-1.1))  and 
(4.201(a) 

shall  hold  and  the  parabolic  problems  (4.1)- (4.4)  has  a  locally  unique  solution  u  which 


ft  is  open,  g^  >  £  >  0  on  R  *  ft 


fulfills 

(4.201(b)  u  e  c3( [0,T] ,  H4( [-1,1] 1)  n  c4( [0,T] ,  tj( (-1,11 )) 

(4.201(c)  u  e  c3(  [--,o] ,  h4(  [-i,m  i,  ux  e  c(  [-«,o] ,  h4(  [-1,1] ) )  o  a  . 

Assumptions  on  the  history  u  and  on  (4.1),  (4.4)  which  guarantee  the  required  smoothness 
of  u  can  be  deduced  from  Renardy  (1981b). 

A  lengthy  calculation  shows  that 

(4.21)  IG.  . (W(N))I  ,  ,  -  0<k2)  +  0(h2) 

k,h  v«X> 


hold.,  where  W(N>  -  ((0^..,)^.,  end 

A-4  (X  -I)2  (X  +1)2 

(4.22)  -  u(xi,tj)  +  - - -  y_(t - - -  y+<tj),  i  -  -M(1)M,j  <  N 

Ai  A1  M  A1 

with  W*:  -  W~_,  ,  W J  *  WJ  .  For  (4.21)  we  used  the  boundary  convergence  result 

M*  1  1  1  “H  *  1 

(3.34)  (convergence  of  the  second  derivative  of  the  boundary  grid  functions  y^, 
yjj  is  necessary  here). 

As  expected  our  scheme  is  consistent  of  second  order  at  the  *  exact'  solution. 
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For  the  stability  analysis  we  calculate  the  Frechet  derivative  of  at  w 

(jet  tiny 

(^.h<s<N>,v<N,,r  -  {+4"vi  -  . 

-  52(xi.tn.h.k)A6v^  -§"3)i(h,k)(^V<N') 


(N) 


( 4*  23 ) (a ) 


( 4*  23 ) (b) 


,(N) ,  Cn 


-^4)i(h.k,(W.v«N,)tn 

for  x  =  -M(1 )M,  0  4  n  <  N 

<G’  «  V?,  i  «  (-M-1 )  ( 1 )  (M+1 ) ,  n  <  0  . 

k,h  ' 


i  i 

Here  V<N>  =  <  6  c’(H2)  such  that  v£+1  -  v"_,  and  V^.,  -  holds 

for  n  <  N.  We  obtain 

N  'n 

(4*24)  g  (x.,t  ,h,k)  =  g  (A6u(x.  ,t  ),(i  Uu(x  ,t.  >>,  >  )  >  e  >  0 

i  1  n  1  in  k  i  j  j  ™ 


*2<xi'tn'h,k>  ”  9ii,A5u(xi'tn)'<ik(4u<xi'tj))j— * 


(4.25) 


N  .  fcn . 


VN  *n 


and  g”  ,  g”  .  are  linear  functionals  on  involving  Frechet  derivatives  of 

(  3)i  ( 4)  1 

g  of  at  most  order  two* 

Similarly  to  the  consistency  result  we  get: 

t  2 

(4*  26)  g.  (x.  ,  t  #h,k)  *  ),u  V/»  +  hV(x  ,t  +  (x,t,h,k) 

l  i  n  1  xt  i  n  x  i  iin  *  i  11 

(4.27)  «2(*i.tnth.k)  «  ♦  h2VVVh,k>  *  Wh'k) 

where  the  vectors  (P#(x. ,t  ,h,k))*  (a# (x, , t  ,h, k ) )*  t  »  1,2  are  uniformly  bounded 

tin  i«*H  tin  i**s 

-2 

in  L  * 
n 

We  investigate  the  linear  equation 

(4.  28)  S  „A“V>  -  r® 

C,  ■  C-  C.  •  C. 

for  n  <  0*  As  in  the  continuous  case  (see  Renardy  (1981b))  we  set 
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(4.29) 


(4.29)  pf  »  Av*1;  Rn  *  G’V';  Qn  »  A  A  V? 

i  xi  x  t  +  -  i 

and  get  the  parabolic  system  of  difference  equations 


(4. 30) (b) 


(4. 30 ) (c) 


(4.30Ma)  «+p"  -  Ar"*1 

(4.30) (b)  «V  -  A  A  Rn+1 

l  +  -  l 

(4. 30 )  (c)  S+r"  ‘T9.lx.it  ,h,k)A  A  (r"+1  +  r")  + 

i  z.  i  i  n  t*  “  l  x 

+  \  *tn#b*K)A(R^+1  ♦  r”)  + 

t  A7J  ......  _  (N)  An  ......  _  ( N )  .Si  n+1 

+  g3(i)(hfk)(xkPi  )  +g4(i)(h,Ic)UKQi  )  ♦  F. 

for  i  »  -M( 1 )M  and  0  <  n  <  N.  For  n  <  0  we  have 


n  i _n  _n  /".n  _n  .  .  n n 

P.  *  AF.  ,  R.  a  0  F.  #  Q,  *  A  A  F 
i  i  x  i  +  -  i 


n  n 
Ru.i 9  *  u 


M+1  M-1 9  -M-1  -M+1 


holds.  Assuming  that  T  is  sufficiently  small  we  regard 


i4.30)(c)  as  perturbation  of  the  time-independent  scheme 


(4.32)(a) 


(4. 32) (b) 


(4. 32) (c ) 
where  *-M-1 


6+r"  *  \  Vy  (x^, 0,0,0)  (A^A^R**1  +  A+A_r”)  ♦ 

+  1  92<Xi'0'0,0)<ARi+1  *  ARi*  +  H¥' '  1  *  0  <  n  <  N 


Ar*  *  Ar1* 


0,  0  4  n  <  N 


Ri  *  1  *  (-M-1  )(1  MM+1  ) 


-0  -0 
R-M4*1  9  RM+1 


R^1  holds.  Denoting  Rn  *  (R^M# . . .  ,R|J, .  •  •  ,R^)  we  write 


(4. 32) (a ), (b) , (c)  in  matrix  form 


B0(h,k)Rn+1  -  B1 (h,k)Rn  +  kHn+1,  R°  »  R° 


where  H 


)  has  been  set.  The  zero-Neumann  boundary  conditions 


(4.32)(b), (c)  are  incorporated  in  the  M  *  M  matrices  BQ(h,k),  Bt(h,k)  (see  (4.51),  (4.52)). 

Proceeding  similarly  to  Varah  (1971a,b)  we  derive 
(4.34)  IB*1 (h,k)l_2  <  Cj 


l(B‘1(h,)t)B1(h,)t))il_2  <  C2,  0  <  i  <  N 


as  h  ♦  0,  k  ♦  0  where  C^,  c2  are  independent  of  h,k.  (4.34)  allows  to  rewrite  (4.33) 


(4.36)  Rn+1  -  C(h,k )Rn  +  kB^Ui.UH"*1,  R°  =  R° 

with  C(h,k)  *  (h.lOB^ (h.k). 

From  (4,35),  (4.36)  and  consideration  similar  to  Benderson  (1971)  we  get  the  stability 


estimate 

(4.37) 


Now  we  assume 
(4.38) 


max  lRnl  +  max  I6+Rnl 

0<n<N  0<n<N  L? 

n  h 


<  c3(Iru 


+  max  I H°l  )  • 

n  0<n*N  Ev 


Xk,  X  =  const  ^  0  as  h,k 


For  the  difference  scheme 


(4. 39)  (a ) 

( 4. 39 )  (b ) 

(4. 39)  (c ) 


and 

( 4.  39 )  (d ) 
(4.39) (e) 


6  P. 
i 


A_n+1 

AR 

i 


6+on  =  A  A  Rn+1 
*i  +  -  i 

6+  n  t  a  ,  „  ^  .  ,  ,„nfl  «.n. 

R^  *  *2  gi  x^»  0, 0,0)  A^A_(R^  +  R^)  + 

1  a  /  ^  n.  n+1 

+  T  92(xi,0'°,0)A(Ri  *  Ri*  *  Hi 


4(1  )M 


0  <  n  <  N 


Ar 


-M 


Arm 


0,  0  <  n  <  N 
% 


-0  -0  *0  -0  *0  -0 

R.  »  R.,  P.  =  P,r  Q.  *  Q,f  i  ~  -H(  1  )M 

i  i  i  i  x  4 


we  get  immediately 

(4.40)  max  lpnl  +  max  2  ♦  max  lgnl  + 

0<n<H  0*n<N  0<n<N 


+  max  I6+Qnl  ♦  max  lRnl  +  max  l5*Rnl  )  < 
0*n<N  0<n*N  0<n<N 

<  const(l?°l  ♦  IQ°I  *  IR°I  +  max  Ih"|  )  . 

S  £h  Hh  0<n<N  S 


We  now  interpret  (4.30)  as  perturbation  of  (4.39).  Therefore  we  write 
(4.41 ) 


3°  (i  p(N))tn 


a0  ..  ACM)  .Si  ^  A0  ..  ^Oo/n 
g3(i)(Vi  >  +93(i)(lkPi  > 


(4.42) 


(N)ven  _  A  0  ,4  a(N)  .  a0  /4  ^(N)  ,fcn 


*4(i)UK«i  ’ 


+^4(i)(ikCi  1 


where  *  0  for  j  <  0;  for  j  >  0  and  ^1  *  *  0  for  3  >  0; 

P^  *  Ri'  ^i  *  ®i  *°r  ^  *  0  and  incorporate  the  second  terms  of  the  right  hand  sides  of 
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(4.4T ),  (4.4 2)  into  the  inhomogeni ty  f"  obtaining  FA  .  The  perturbation  originating 
from  taking  ^4(i)  instedd  of  ^3(i)  and  ^4(i>  re8p>  wiU  66  investigated 

later.  We  denote  the  sum  of  the  two  remaining  terms  in  (4.41),  (4.42)  by  (P-1 , ) "=1 
and  from  (4. 20) (a)  we  derive  the  estimate 

(4.43)  MJMpJ.O?  >?,,>?,  J  2  <  const<  max  1^1  ,  +  max  Iff*!.) 

i  i  i  i— M  l2  i«j<n  Lf  1<j<n  LT 

h  h  h 


uniformly  as  h,k  ♦  0. 

We  denote 
(4.44) 


<*>.„ . ^m'c-m . cm'r"« . rm> 


and  rewrite  the  perturbed  system  (4.30)  (with  tR  ■  0  in  g^,  g2  end  g^^*  ^4(i) 


instead  of  g,f . g 


(4.45)  (a) 

(4.45 )  (b) 
where 


3(i  ) '  ^4(i  ) 


)  as 
.  n+1 


E0(h,)t)An+1  »  E1  (h,k)An  +  XLQ(h,k)An+1  +  kL^h.kJA"  ♦  kOtA3)^  +  kF1 


n+1 


(4.46) (a) 


.(h.k) 


2M+1 


2M+1 


-kD^h.k) 

-kD2(h,k) 


}  2M+1 

}  2M+1 

}  2M+1 


( 4. 46) (b) 


Ej  (h,k) 


6 

e 

(h,k) 


( 4.46) (c ) 


a<Aj>;a1 


(d)  F 


n-M 


0 
0 

(F"+V4  „ 

i  i*-M 
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holds*  Here  we  denoted  the  difference  rad trices 


<4.47)<a) 


( 4. 47 ) (b ) 


The  first  and  last 
Neumann  conditions 
The  matrices 
instead  of 
(4.40) 


Vh'k)  *  2h 


10-1  0 

\  \  \ 

1  0  -1 

\  \  \ 

e  10-1 

ooo 


2K+1 


D,(h,k)  *  — r 
4  fl 


-2  2  0 
1-219 

\  \  \ 

1  -2  1 

W  \ 


1  -2  1 

0  2  -2  J 


2M+1 


rows  of  D,(h,lc),  D2(h,k)  come  from  incorporating  the  homogeneous 

L0(h,M,  LjOijk)  are  derived  by  taking  g1  (xt ,0, h,k ) ,  g2(xif0,h,k) 
►,  0, 0 )  and  g^ix^, 0,0,0)  reap •  Because  of  (4.26)  ,  (4.27)  we  get 
lLA(h,k)l_2  <  const.  ,  %  -  0,1 


uniformly  as  h,k  ♦  0, 


We  easily  obtain 


(4,49) 


E~’(h,lc) 


I  8  kDj (h,k)B~' (h.k) 
8  I  kD^Oi.kjB”1 (h.k) 
8  8  B~' (h.k) 


and 


(4.50) 


for  i  >  0. 


(E'V.kJB,  (h.k))1 


X 

I  8  -kD  (h.k)  l  C(h,k)3 
j-1 
i 

0  I  -JcD2(h,k)  l  C(h?k) 
j-1 

0  0  C(h,k)L 


A  3 

We  set  -  gt  (xi#0f  0#0)#  0^  -  ^  (xi  ,0, 0,0 )  and  get  from  (4.32)  (a),  (b) 


(4.51) 

(4.52) 


B0(h,k)  -  I-(|  A2D2(h#X)  +|A1D1(hrk)) 
k  *  _  .  .  k 


(h,k)  -  z+f  A2D2(h'k>  +^A1ol(h,k)) 

where  A2  -  diagfa^,. . . ,oQ,.  •  •  /aM )  and  A?  -  diagtfl^, . . .  ,30».  • .  ,0^)  holds.  Since 
2 


h~  »  Xk  we  get 

(4.53) 

Obviously 

(4.54) 

solves  the  matrix  iteration 

(4.55)  (a) 


lE  ^hjk)!  <  const. 
Yl  -  k  l  C(h.k)1 

j-1 


(4.55)  (b) 

From  (4.37)  we  get 

(4.56) 


Y1*1  -  C(h,k)Yi  +  kC(h.k)  ,  0  <  i  <  N-1 

Y°  -  8  . 


max  It  I  <  const. IB  (h,k)l  • 
0<i<N  i* 


.9 

Because  of  (4.38)  B^(h.k)  is  uniformly  bounded  (in  L^)  and  therefore  the  right 

hand  side  of  (4.56)  is  uniformly  bounded  as  h,k  ♦  0.  We  obtain 
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(4,57) 


'V  \  Ci'r2' 


A. 

ID  Jc  J  C"*  I  <  const., 
3=1  ^ 


0  4  l  <  N 


uniformly  as  h,k  ♦  0. 

It  is  interesting  to  note  that  the  mesh-si2e  restriction  (4.38)  is  not  necessary  for 
the  implicit  Euler-scherae  since  for  this  fully  implicit  scheme  8^  »  I  holds* 

(4*57)  implies  L2-stability  of  (4.45): 

14.58)  I  (E^hfkjE  (h,*)))1!  4  const,  for  0  4  i  4  n 

<lV 

h 


uniformly  as  h,k  ♦  0. 

By  proceeding  similarly  to  Richtmeyer  and  Morton  0985,  Chapter  3.9)  we  get 
(4.59)  IAnl  ,  ,  <  const.  (IA°I  ,  .  +  max  l(F?)*  I  ,) 

‘V  (V  1<:i<n  5. 

where  (An)N  _  solves  (4.45).  For  the  (first)  time  difference  quotient  of  Rn  we  get 
n*0 


S*Rn  -  ±  (C(h,k)  -  l)Rn  +  (Bp1(h,k)L0(h,)t)Rn  + 


(4.60) 


♦  B-1<hflc)I1<h.)c)Rn+1  +  8o1<h',tHai(I>i'ei)j  =  1,i=-M  + 
where  L^(h,k),  i  *  0,1  stands  for  the  (block)  matrix  in  the  (3.3)  position  of  La(h,k). 
From  (4.36),  (4.37)  we  get 

(4*61)  lfi+Rnl  4  const.  (IR°I  +  max  lA^I  +  max  I (F? )M_  I  ) 

Ll  Hi  1<j<n  (L2)3  i<j<n-f  1  1  1=~M  L2 

h  h  h  n 


Similarly  we  get  bounds  for  the  spatial  difference  quotients  and  using  (4.59)  we  obtain  the 
stability  estimate 
(4*62) 


lRnl  <  const.  <IR°I  *  IP°I  ,  +  IS°I  _  +  max 

„2  „2 


^  1  i=_M  ^ 


From  (4.29)  we  get  immediately 
(4*63) 


KG'  (Vi(N)))_,l  ,  ,  ,  <  const. 

Ck<Hh*  n*  *  Ck(V 


uniformly  as  h,k  *  0  where 
(4.64) 


«,N>  - 
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t 

i, 


and 


(4*65 )( a) 
( 4.65 ) ( b) 
holds. 


Wi  =  '  1  *  (“M*1 ) ( t ) (M+1 ) ,  j  <  0 

i  -  (-M-IKtHM+1),  0  <  j  <  N 


=  W° 


As  a  lengthy  (but  easy)  calculation  shows  G^  ^  is  uniformly  Lipschitz  continuous  in 

a  sphere  S^tW^)  <=■  C^(H^)  whose  radius  P  is  independent  of  h,k: 

/>.  v  ,„(N).  ,„(N)._  *  iv(N)  _  ( N ) « 

(4.66)  IG*  (Y  )  -  G'  (Z  )l  .  0  0  -  <  const. lY  -  Z  •  , 

k'h  k'h  cW>*C.(H2.i2>  c’(H2) 

k  n  k  n  n  k  n 


for  Z(N),  Y(N)  e  Sp(W(N)).  Since  (4.20)  implies  that 


(4.67) 


.w(N)  -a(N)« 


o(l)  as  T  ♦  0 


holds,  the  estimate  (4.63)  is  also  fulfilled  by  G*  (W^N^)  if  T  is  sufficiently  small. 

k  ih 

Applying  Keller*s  (1975)  theory  we  obtain 

Theorem  4.1.  Assume  that  the  assumptions  (4.20)  and  (4.38)  hold.  Then  for  h,k 


(4.68) 


sufficiently  small  the  scheme  (4.8),  (4.9),  (4.10),  (4.11)  has  a  locally 
unique  solution  ((U^)^*w„|  if  T  is  sufficiently  small  and  the 

convergence  estimate 

'i-m'v  ,  “0(h2> 

ic  n 


as  h  ♦  0  holds.  The  (abstract)  Newton  method  for  (4.13)  converges  quadratically  from  a 
(sufficiently  small)  sphere  of  starting  values  whose  radius  is  constant  as  h,k  ♦  0, 

Of  course,  T  can  be  taken  independently  of  h,k. 

Remark .  Renardy  (1981)  assumed  that  g  s  R  *  fl  ♦  H*( 1-1,1]),  c  H^([-1,1])  holds 
instead  of  (4.20) (a).  This  is  a  more  realistic  assumption  (with  respect  to 
viscoelastic  problems)  since  H1((-1,U)  is  a  Banach  algebra  (elements  can  be 
multiplied),  but  the  perturbation  approach  (4.45)  would  not  go  through  as 
presented.  However,  this  is  a  technicality  which  can  be  repaired  by 
incorporating  one  more  x-difference  quotient  into  the  spaces. 
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As  an  example,  we  apply  our  Crank-Nico lson-type  scheme  to  a  model  for  the  stretching 
of  a  thin  filament  of  a  viscoelastic  fluid  when  a  force  f  is  applied  to  its  ends  (derived 
by  Renardy  (1961b)): 

(4.69)(a)  Pu  *  (3n  +  u  Ftu1')  -  GfuS),  x  6  [-1,1],  t  6  [0,T] 

tt  *  DC  X  2  X 


( 4*69) (b ) 

(4.69)  (c) 

(4.69 )  (d ) 
where  y(t) 


u(x,t)  *  u(x,t) 
ux(1,t)  *  y(t) 
u^t-1 ,t)  *  y(t) 
solves 


x  €  (-1,1],  t  €  {-*>,0] 
t  G  [ 

t  6  [-«,T] 


(4.70) (a ) 


31  +  y(t)F(yt)  -  — - —  G(yt) 

y  (t)  y*<t> 


f(t),  t  6  [0,T] 


(4.70) (b)  y (t )  -  u  (1,t)  *  u  (-1,t)  *  t  6  (-*,0)  . 

The  history  u(x,t)  is  assumed  to  be  an  odd  function  of  x  for  all  t  6  (-*,0),  F,  G 
are  functionals  (of  Volterra  type).  The  parameters  P,  n  and  the  physical  meaning  of  u 
is  explained  in  Renardy  (1981b). 

At  first  we  apply  the  midpoint  rule  (3.3),  (3.4)  to  the  boundary  problems  (4.69)  and 
obtain 


n  n+1 


( 4.71 ) (a)  12h  -  n  -  F((i  y(N))  n*1/2)  - 

(y  +Y  ) 


7  n'Vvra  G(liky(N,)  n+V2>  -  f<t„+1/2>'  0  <  n  <  ■ 

(y  +y  ) 


(4.71  Mb)  yn  ■  u  <1,t  ),  n  <  0  . 

x  n 

(4. 8) -(4.11)  applied  to  (4.69)  gives 


(4.72) (a) 


«A  A  u"(Au")2  -  2A«U?Au"a  A  u” 

.  3n  — - - i - - — i-Z-Z  i. 

(At)")4 


t 
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| 


♦  Au^  n)(i)c\A.UiN,>  "  + 


Al^A  A  Un  t 

♦  2  —  — t  ~-  *-  G((iAu[N))  ")  - 

<Au">4  k  1 


^v«v<>>tn,‘w-“r>tn 


for  i  <  0  <  n  <  N 

(4.7 2) (b)  U?  -  u(x.  ,t  ),  i  •  n  <  0 

i  in 

(4.72)  (c)  Au£  -  yn  ,  n  <  N 

(4.72) (d)  Al/^  «  yn  ,  n  <  N  . 

Here  F-| ,  G1  denote  the  first  Frechet  derivatives  of  F  and  G  reap. 
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